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covariant gauge using a mixture of analytical and numerical methods. We improve the precision of 
previous results. 
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I. INTRODUCTION 

Muonium is the hydrogenlike bound state of a positive 
muon and an electron. Unlike hydrogen, or any other 
bound state involving hadrons, muonium is free from the 
complications introduced by the finite size or the inter- 
nal structure of any of its constituents. Therefore, it al- 
lows for a very precise test of bound-state QED, and can 
be used to restrict models of physics beyond the Stan- 
dard Model. Measurements of the ground-state hyper- 
fine splitting of muonium are used to extract the muon 
to electron mass ratio m^/rrie and the muon to proton 
magnetic moment ratio ^i^J ^p |1!|. The value of ^i^J ^p 
is required for obtaining the muon anomalous magnetic 
moment from experiment [2]. In addition, the hyperfine 
splitting can also be used to determine the fine structure 
constant a. For a review of the present status and recent 
developments in the theory of light hydrogenic atoms, 
see [3|,||. 

The leading-order hyperfine splitting is given by the 
Fermi energy Ep (defined in Eq. ([T])). Its correc- 
tions are organized as a perturbative expansion in pow- 
ers of three parameters: Za, describing effects due to 
the binding of an electron to a nucleus of atomic num- 
ber Z\ a (frequently accompanied by I/tt) from electron 
and photon self-interactions; and the ratio of electron 
to nucleus masses, m/M . The main theoretical uncer- 
tainty comes from three types of yet unknown correc- 
tions: single-logarithmic and non-logarithmic corrections 
of order a{Za)'^{m/M)Ep, and non-logarithmic correc- 
tions of order a'^{Za){m/M)EF and {Za)^{m/M)EF 
(some terms are known for the first case [51 ). 

In this paper we focus on the second-order radiative- 
nonrecoil corrections to the hyperfine splitting (of order 
a'^{Za)Ep). The total result for these corrections was 
found by Eides and Shelyuto [Q] and Kinoshita and Nio 
Q- Our result improves their precision by over an order 
of magnitude. Our central value is slightly lower than, 
but compatible with, that of p. 

In Sec. |TT]we present the details of our approach, and 



in Sec. IIIII we present our results. In Appendix |^ we 
show analytic results for two sets of diagrams. 



II. EVALUATION 

We consider an electron of mass m orbiting a nucleus of 
mass M and atomic number Z . In this paper we consider 
the nucleus to be a muon, but we will keep Z explicit in 
order to distinguish between the binding contributions 
{Za) and the radiative ones (a). 

We are interested in corrections to the hyperfine split- 
ting of the ground state of muonium of order a'^{Za)Ep 
and leading order in m/M , where 
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Here g is the gyromagnetic factor of the nucleus^ (in our 
case, a muon, but our final result in Eq. (|19p applies 
to any hydrogenlike atom). In order to compute these 
corrections, we consider the scattering amplitude 



iM = [u{p)Qiu{p)][v{P)Q2v{P)] 



(2) 



where u[p) is the spinor for the electron, v{P) is the 
spinor for the muon, p = (m, 0) and P = (M, 0). Qi and 
Q2 are given by the Feynman rules describing the sum 
of the diagrams shown in Figs. [1] and [2] In these figures, 
the sum of the direct and crossed interactions between 
the electron and the muon is represented by a dotted 
line, as shown in Fig. [3l We define a bound-state wave 
function ^ = uv, so that Eq. ^ becomes 



iM = -Ty{^^ QiijQ2} ■ 



(3) 



^ It includes the corrections from the anomalous magnetic mo- 
ment, which factorize with respect to the corrections considered 
in this paper. This is no longer true when considering non-recoil 
corrections. See e.g. Q. 
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Figure 1: The different sets of vacuum polarization diagrams 
and light-by-light diagrams (set IV). Each set represents the 
drawn diagram plus all the possible permutations of its pieces. 



Depending on the relative alignment of the spins of the 
constituent particles, an S state can either belong to the 
J = 1 triplet or the J = singlet. The triplet and singlet 
states are often denoted by the prefixes ortho- and para-, 
respectively, and their wave functions are given by [8| 
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where ^ is the polarization vector. We average over the 
directions of £, by considering the four- vector ^ = (0,^) 
and using the identity 
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We use dimensional regularization with d = 4— 2e dimen- 
sions. Thus, an important issue is the definition of 75, 
which is an intrinsically four-dimensional object. Since 
we do not have to evaluate traces with an odd number of 
75 matrices, we can treat them as anticommuting. 

The energy shift created by the radiative corrections 
depicted in Figs. [T] and [2l for either the singlet or triplet 
configurations, is given by 
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where |-0n(O)P — (Zajj.)'^ /{Trn^) is the squared modulus 
of the wave function of a bound S state with principal 
quantum number n and reduced mass fi. The hyperfine 
splitting is then simply 



SEhfs = SE. 
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In order to evaluate the loop integrals represented by 
the Feynman diagrams we use the method of regions 9] to 
construct an expansion in the small ratio m/M . There 
are several possible contributing regions, where one or 
more of the loop momenta scale like m or M . However, 
we are only interested in the leading order in m/M , which 
is given by the region where all loop momenta scale like 
m. If fc ~ m, we can expand the contribution from the 
muon line in the sum of the direct and crossed diagrams 
of Fig. 1 
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Figure 2: The diagrams involving a two-loop electron self- 
interaction and vertex corrections. 
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Figure 3: The sum of the direct and crossed diagrams is rep- 
resented by a dotted line (the double line represents the prop- 
agator of the muon) . 



where 
Ti = 

T2 = 
^3 = 
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We used the equation of motion to set some terms in 
the numerator to zero, and we arranged the terms in the 
expansion in such a way that the three different Dirac 
structures that are important for the calculation of the 
hyperfine splitting appear explicitly. We will now see 
that only T2 can contribute to the splitting. 

Consider the Dirac structure of ipTi and anticommute 
the gamma matrices, for both para and ortho states: 
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Now we can write 



ex Tr{V'tQiXTj=Tr{xTi^^Qi}. (15) 



Using the expressions in Eqs. p^ and (jTl)) it is easy 

ort/io (after averaging 



to see that x^^^V'^a™ = X??;*'"'V',^ 



over polarizations). This means that Ti gives the same 
contribution for para and ortho states. Therefore, when 
we subtract these contributions in order to compute the 
hyperfine sphtting, they cancel out. 

If we consider T2 instead, defining XT2 in analogy 
with Eqs. p^ and (ITU) , we can see that K^^^'^i^para 7^ 

X'T*'^°'Vortho^ ^^ ^^i^ term will not cancel in the subtrac- 
tion. The difference between the para and ortho states 
comes solely from terms in Xt^^" ^^^t are totally anti- 
symmetric in a and /3. Therefore, when we consider the 
Dirac structure of T3, which is but a symmetrization of 
that of T2, these terms will vanish, and so T3 will give no 
contribution to the hyperfine splitting either. 

Thus, we have seen that the only term that contributes 
to the hyperfine splitting is 
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This is valid at all orders of alpha, and all orders in m/M . 
We can then substitute the scalar part of the nucleon 
propagator by a Dirac delta in all our calculations, since 
we are only interested in the leading order in m/M and 
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We used dimensional regularization, and renormalized 
our results using the on-shell renormalization scheme. 
For all the photon propagators in Figs. [1] and [H we used 
a general covariant R^ gauge. The overall cancellation 
of the dependence on the gauge parameter in the final 
result provides us with a good check for our calculations. 

We used the program qgraf 10] to generate all of the 
diagrams, and the packages q2e and exp jlll . Il2| to ex- 
press them as a series of vertices and propagators that can 
be read by the FORM ^13,] package MAT AD 3 ^^. Finally, 
MATAD 3 was used to represent the diagrams in terms of a 
set of scalar integrals using custom-made routines. In this 
way, we represented the amplitude Ai in terms of several 
thousand different scalar integrals. These integrals can 
be expressed in terms of a few master i nteg rals by means 
of integration- by-parts (IBP) identities [3. We used the 
so-called Laporta algorithm |16l. |l7ll as implemented in 
the Mathematica package FIRE [lg|, to reduce the prob- 
lem to 32 master integrals. The master integrals for this 
calculation are the same ones we found in [l^. All defi- 
nitions and results for the integrals can be found in this 



Table I: Comparison between our results for sets of diagrams 
of Fig. [T] and those of [23 - [24| . Numbers ending in an ellipsis 
indicate an analytic result, which we show in Appendix \K\ 
No error was given for the numerical result of set I in [23| . 



Set 


This paper 




Refs. [22-241 


I 


-0.3107420427660S 


(3) 


-0.310742 


II 


-0.668915 . . . 




-0.668915 . . . 


III 


1.867852 . . . 




1.867852 . . . 


IV 


-0.4725146(2) 




-0.472514(1) 


V 


36/35 




36/35 



reference. However, one change was made for this calcu- 
lation. In order to obtain better numerical precision, we 
performed a change of basis, so that instead of working 
with /i4 = F(l, 0, 0, 0, 1, 1, 1, 1) we worked with 

^(1,0,0,0,1,1,1,2) 
= 44.55822275(2) - 427.382296(2)e + 0{e^) , (18) 

which was obtained using the Mathematica package 
FIESTA 1.2.1 [20] with integrators from the CUBA li- 
brary 2l|. 



III. RESULTS 



Our final result for the hyperfine splitting is 



SEhfs = 0.77099(2) 
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This correction was also found by Fides and Shelyuto 
and by Kinoshita and Nio (7|. Their results are 
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Our result is a little over one order of magnitude more 
precise than that of [yj , and almost three orders of mag- 
nitude more precise than the one in [7| . Our central value 
is slightly lower than in Q, by about 1.5cr (taking as a 
the larger error). It agrees with the result of Q within 
its much larger error estimate. For the ground state of 
muonium, our result reads 
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0.42524(1) kHz. 
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We compared our results for the individual diagrams 
and those found in the literature [0, 0, I22H24J . Our re- 
sults for the gauge-invariant sets of diagrams of Fig.[T]are 
presented in Table ID For the diagrams of Fig. [5] we chose 
the Fried- Yennie gauge [23, l26| , in which all diagrams are 
infrared finite. Our results are presented in Table HIl 

The sum of all central values in the second column of 
Tables U and ini gives the coefficient 0.77099 in Eq. (E 



Table II: Comparison between our results for diagrams as 
(in Fried- Yennie gauge) and those of [g]. 



Diagram 


This paper 


Ref. f6] 


a 


9/4 


9/4 


b 


-6.6602948853575169751(3) 


-6.65997(1) 


c 


3.9324055550472089860(4) 


3.93208(1) 


d 


-3.9032816968990(2) 


-3.903368(79) 


e 


4.5667195410288(2) 


4.566710(24) 


/ 


-37rV8 + 19/64 


-3.404163(22) 


9 


7rV2 - 9/4 


2.684706(26) 


h 


33/16 


33/16 


i 


0.05454(1) 


0.054645(46) 


J 


-7.14963(2) 


-7.14937(16) 


k 


1.4658690989997(5) 


1.465834(20) 


I 


-1.98334(3) 


-1.983298(95) 


m 


3.16949(2) 


3.16956(16) 


n 


-3.59661163(2) 


-3.59566(14) 





1.80476(5) 


1.804775(46) 


P 


3.507035(6) 


3.50608(16) 


q 


-0.80380(3) 


-0.80380(15) 


r 


1.05247(3) 


1.05298(18) 


s 


0.277336777308(2) 


0.277203(27) 



The reason for the discrepancies (and their cancella- 
tions) is most likely the different treatment of infrared 
divergences in Q and this paper. In Q, the Fried- Yennie 
gauge was set from the beginning, and all spurious in- 
frared divergences were canceled before the integration 
over the diagram's loop momenta, which was performed 
in four dimensions. In our calculation, we used a gen- 
eral gauge parameter, and dimensional regularization to 
deal with infrared divergences, which would only vanish 
after setting the gauge in the final expression. As noted 
in [23], there is a difference between setting the Fried- 
Yennie gauge and sending the infrared regulator to zero 
before or after integration. It is not surprising then that 
we obtained different results than [0] for gauge-dependent 
diagrams, but that most of the differences cancel in the 
final, gauge-invariant result, making it compatible with 
the previous calculation. 

Using the setup of the calculation of the hyperfine split- 
ting one can also find the Lamb shift, as it is given by 

CJ-, SEortho {d— I) + SEpara ,^„-. 

0-CjLamb = "1 ■ [^o) 

We obtained in this way the same results as in [19,]^. 
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The error of that result is however not obtained from the 
sum of the errors of the diagrams in the tables. Once 
we decompose the problem into the calculation of mas- 
ter integrals, the diagrams arc no longer independent, as 
the same master integral contributes to several different 
diagrams. Thus, to find the error of our total result, we 
first sum all diagrams and then sum all the errors of the 
integrals in quadrature. 

We found new analytic results for diagrams g and /, 
shown in Table |TT1 For completeness, the known analytic 
results for sets II and III of the vacuum polarization di- 
agrams are given in Appendix |^ as well. 
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Appendix A: Analytic results 

Here we show the analytic results for sets II and HI of 
the vacuum polarization diagrams, found in [2^, 
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Set HI 
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—23/378. This does not affect any of the other results presented 
in that paper. 
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